#ifndef TOPOLOGY_HPP
#define TOPOLOGY_HPP

#include "meshTypes_fwd.hpp"
#include "mesh.hpp"
#include <array>
#include <vector>

namespace hms
{

namespace aux {
	Index otherEntry(
		const Ref<const Array2i>& col,
		Index entry
	);
}

/* if you have an edge consisting of vertices A and B, then specifying A as the
 * current vertex returns B, while B returns A. */
template<typename EdgeType>
Index otherVertex(
	const EdgeBase<EdgeType>& edge,
	Index currentVertex
);

template<typename MeshType>
Index otherVertex(
	const MeshBase<MeshType>&,
	Index edge,
	Index currentVertex );

Index otherVertex(
	const Array2Xi& edgeVertices,
	Index edge,
	Index currentVertex );



/* takes two edges, returns the index of the shared/connecting vertex */
Index sharedVertex( const Array2Xi& edgeVertices, Index edgeA, Index edgeB );
Index sharedVertex( const Mesh&, Index, Index );
Index sharedVertex( const Edge&, const Edge& );

template<typename EdgeType>
Index sharedVertex( const EdgeBase<EdgeType>&, const EdgeBase<EdgeType>& );



/* returns all edges connected to each vertex */
std::vector<IdxVector>  edgesAtVertices(
	Index nVertices, const Array2Xi& edgeVertices
);

std::vector<IdxVector>  edgesAtVertices( const Mesh& );



std::vector<IdxVector> findLoops(
	const Array2Xi& edgeVertices,
	Index edge,
	const std::vector<IdxVector>& edgesAtVertex );

std::vector<IdxVector> findLoops(
	const Mesh& mesh,
	Index edge,
	const std::vector<IdxVector>& edgesAtVertex );

std::vector<IdxVector>  findLoops(
	const Mesh&,
	const Edge&,
	const std::vector<IdxVector>& edgesAtVertex );



ArrayXi edgeLoopToVertices(
	const Array2Xi& edgeVertices,
	const Ref<const ArrayXi>& edgeLoop
);

ArrayXi edgeLoopToVertices(
	const Mesh&,
	const Ref<const ArrayXi>& edgeLoop
);

/* takes an edge loop and returns the unique vertices */
ArrayXi edgeLoopToVertices( const Mesh&, const IdxVector& edgeLoop );

ArrayXi edgeLoopToVertices( const Mesh&, Index cell );


enum class Side
{
	left,
	right,
	on,
};

Side whichSide(
	const Ref<const Vector2s>& a,
	const Ref<const Vector2s>& b,
	const Ref<const Vector2s>& p );

Side whichSide( const Edge&, const Ref<const Vector2s>& );
Side whichSide( const Mesh&, Index, const Ref<const Vector2s>& );

Index findBottomLeft( const Ref<const Matrix2Xs>& vertices );

Index findBottomLeft(
	const Mesh&,
	const ArrayXi& vertices
);


/* vertices may be ordered clockwise or counter-clockwise. Both start with the
 * bottom left vertex. */
enum class VertexOrdering
{
	clockwise,
	counterclockwise,
};

void orderVertices(
	const Matrix2Xs& vertexCoords,
	Ref<ArrayXi> vertices,
	VertexOrdering=VertexOrdering::counterclockwise );

void orderVertices(
	const Mesh& mesh,
	Ref<ArrayXi> vertices,
	VertexOrdering=VertexOrdering::counterclockwise );



void orderEdgeLoop(
	const Mesh& mesh,
	const Ref<const ArrayXi>& orderedVertices,
	IdxVector& edgeLoop );

void orderEdgeLoop(
	const Mesh& mesh,
	const Ref<const ArrayXi>& orderedVertices,
	Ref<ArrayXi> edgeLoop );



/* takes a list of vertices and returns the interconnecting edges */
IdxVector connectingEdges(
	const Mesh&,
	const IdxVector& vertices,
	const std::vector<IdxVector>& edgesAtVertex );


/* returns the bounding box of a polygon as bottom left and top right corner */
AlignedBox boundingBox( const Ref<const Matrix2Xs>& );
AlignedBox boundingBox( const std::vector<Vector2s>& );

/* Determines whether a point lies in a polygon. Uses winding number algorithm,
 * thus works for simple and complex polygons. Points on edges are not included.
 * For repeated testing of the same polygon, it is recommended to check the
 * bounding box first. */
bool inPolygon( const Ref<const Matrix2Xs>&, const Ref<const Vector2s>& );
bool inPolygon( const std::vector<Vector2s>&, const Ref<const Vector2s>& );

/* takes a list of vertices as a polygon and returns the vertices of the
 * passed-in mesh which lie inside that polygon */
IdxVector verticesInPolygon( 
	const Mesh&,
	const std::vector<Vector2s>& polygon );

/* takes a list of vertices forming a non-intersecting polygon and returns 
 * the edges of the passed-in mesh which lie fully inside the polygon.
 * Calls edgesAtVertices, if last parameter is empty. For repetitive calls,
 * call edgesAtVertices first and pass it in */
IdxVector edgesInPolygon( 
	const Mesh&,
	const std::vector<Vector2s>& polygon,
	const std::vector<IdxVector>& edgesAtVertex = {} );

/* returns cells whose cell centroid lies inside the specified polygon */
IdxVector cellsInPolygon( const Mesh&, const std::vector<Vector2s>&);

/* creates an axis-aligned rectangle from two diagonally opposed corners a and b
 * and returns a list of the vertices of the passed-in mesh which lie inside */
IdxVector verticesInBox(
	const Mesh&,
	const Vector2s& a,
	const Vector2s& b );

/* creates an axis-aligned rectangle from two diagonally opposed corners a and b
 * and returns a list of the edges of the passed-in mesh which lie inside */
IdxVector edgesInBox(
	const Mesh&,
	const Vector2s& a,
	const Vector2s& b,
	const std::vector<IdxVector>& edgesAtVertex = {} );

/* creates an axis-aligned rectangle from two diagonally opposed corners a and b
 * and returns a list of the cells of the passed-in mesh whose centroids lie
 * inside */
IdxVector cellsInBox(
	const Mesh&,
	const Vector2s& a,
	const Vector2s& b );

/* traces the specified vertices and collects all the edges that lie between 
 * them. If they are ordered, no duplicates occur. GIGO applies */
IdxVector verticesToEdgeLine(
	const Mesh& mesh,
	const IdxVector& vertices, 
	const std::vector<IdxVector>& edgesAtVertex );

template<typename MeshType>
std::optional<Index> containingCell( const MeshBase<MeshType>&, Array2s coords );

template<typename MeshType>
IdxVector cellsAlongLine( const MeshBase<MeshType>&, Array2s beg, Array2s end );

} // end namespace hms

#include "topology.tpp"

#endif
